%% Forecast ...
% by IMF Institute for Capacity Development - Course on Monetary and Exchange Rate Policies

%% Clear the workspace
close all;
clear all;
disp('Working ...');

%% Read the model
readmodel;

% Loads the historical data from the file 'history.csv'
% (the file 'history.mat' is saved by the 
% 'makedata' program) into the object 'h'
h = dbload('history.csv');

% Defines the time frame of the forecast
% Don't forget to change the time frame depending on you data and forecast
% period
startfcast = qq(2013,2);
endfcast = qq(2022,4);
fcastrange = startfcast:endfcast;

% Shock to risk premium
% h.shk_prem = tseries(startfcast:startfcast+1,5);

% Shock to aggregate supply (e.g., increase in indirect taxes)
% h.shk_dot_cpi = tseries(startfcast,1); 

% Shock to aggregate demand (e.g., sharp cyclical slowdown)
% h.shk_lgdp_gap = tseries(startfcast:startfcast+3,-0.75); 

% h.shk_lx_gdp_gap = tseries(startfcast:startfcast,-1);

% Command 'simulate' simulates the model 'm' based on the database 'h' over
% the forecast range 'fcastrange'. Results are written into the object 's'
s = simulate(m,h,fcastrange,'anticipate=',false,'nonlinearise=',30,'maxiter=',1000);

%% In the forecast simulation variables can be set as numbers (for instance
% for external development as opposed to automatically used AR(1) process). 
% These preset variables are in the historical database and they require 
% creation of a 'plan' for this variable. The command 'plan' creates plan 
% for model 'm' over the forecast horizon in the object 'simplan'. These
% variables must be marked using the command 'exogenize', while the residual
% used for equalizing the model equation for the fixed variable must be marked using
% command 'endogenize'. Note that any variable that is fixed must have residual term
% as a part of its equation. Finally, the forecast is simulated in a roughly 
% same way as without the fixing. See example below where external inflation 
% and external interest rates are taken, say, from IMF's WEO or Consensus
% Forecast:

% %Defining foreign inflation
% h.dot_x_cpi(startfcast:startfcast+1) = 4;
% % % 
% % % Fixing foreign inflation 
% simplan = plan(m,fcastrange);
% simplan = exogenize(simplan,{'dot_x_cpi'},startfcast:startfcast+1);
% simplan = endogenize(simplan,{'shk_dot_x_cpi'},startfcast:startfcast+1);
% s = simulate(m,h,fcastrange,'plan',simplan,'nonlinearise=',30,'maxiter=',1000);

% %Defining domestic policy rate
% h.rn(startfcast:startfcast+1) = 0.5;
 
%Fixing the policy rate for one quarter
%  simplan = plan(m,fcastrange);
%  simplan = exogenize(simplan,{'rn'},startfcast:startfcast);
%  simplan = endogenize(simplan,{'shk_rn'},startfcast:startfcast);
%  s = simulate(m,h,fcastrange,'plan',simplan,'nonlinearise=',30,'maxiter=',1000);

% % Defining foreign policy rate
% h.x_rn(startfcast:startfcast+3) = 3;
%  
% %Fixing foreign interest rates
% simplan = plan(m,fcastrange);
% simplan = exogenize(simplan,{'x_rn'},startfcast:startfcast+3);
% simplan = endogenize(simplan,{'shk_x_rn'},startfcast:startfcast+3);
% s = simulate(m,h,fcastrange,'plan',simplan,'nonlinearise=',30,'maxiter=',1000);

% % Defining foreign output gap
% h.lx_gdp_gap(startfcast:startfcast+1) = -1;
%  
% %Fixing foreign output gap
% simplan = plan(m,fcastrange);
% simplan = exogenize(simplan,{'lx_gdp_gap'},startfcast:startfcast+1);
% simplan = endogenize(simplan,{'shk_lx_gdp_gap'},startfcast:startfcast+1);
% s = simulate(m,h,fcastrange,'plan=',simplan,'nonlinearise=',30,'maxiter=',1000);

% Command 'dbextend' puts together the historical database 'h' and the
% results of the simulation saved in object 's'. Single database 's' is
% created
s.rn = s.rs;
s = dbextend(h,s);

% Results are saved in file 'fcastdata.csv'  
dbsave('fcastdata.csv',s);

%% Report
% Prepares the forecast report: graphs and tables in the pdf format
Tablerng = startfcast-3:startfcast+7;
Plotrng = startfcast-3:startfcast+11;
Histrng = startfcast-3:startfcast-1;

x = report.new('Georgia - Forecast','visible',true);

% Figures
sty = struct();
sty.line.linewidth = 1.5;
sty.line.linestyle = {'-';'--';':'};
sty.line.color = {'k';'k';'r'};
sty.axes.box = 'on';
sty.legend.location = 'southOutside';
sty.legend.orientation = 'horizontal';

x.figure('Forecats - Main Indicators','subplot',[3,2],'style',sty,...
  'range',Plotrng,'dateformat','YYYY:P');

x.graph('Inflation','legend',true);
x.series({'q-o-q','y-o-y','Target'},[s.dot_cpi s.dot4_cpi s.target]);
x.highlight('',Histrng);
  
x.graph('Nominal Interest Rates','legend',false);
x.series('',[s.rn]);
x.highlight('',Histrng);
    
x.graph('Nominal Exchange Rate - GEL/USD','legend',false);
x.series('',[exp(s.ls/100)]);
x.highlight('',Histrng);
  
x.graph('Nominal Exchange Rate (q-o-q)','legend',false);
x.series('',[s.dot_s]);
x.highlight('',Histrng);

x.graph('Output Gap','legend',false);
x.series('',[s.lgdp_gap]);
x.highlight('',Histrng);
    
x.graph('Real Interest Rate and Real Exchange Rate Gap','legend',true);
x.series({'RIR Gap','RER Gap'},[s.rr_gap s.lz_gap]);
x.highlight('',Histrng);

x.pagebreak();

% Tables
TableOptions = {'range',Tablerng,'vline',startfcast-1,'decimal',1,'dateformat','YYYY:P',...
    'long',true,'longfoot','---continued','longfootposition','right'};

x.table('Forecast - Main Indicators',TableOptions{:});

x.subheading('');
  x.series('CPI ',s.dot4_cpi,'units','% (y-o-y)');
  x.series('',s.dot_cpi,'units','% (q-o-q)');
  x.series('Target',s.target,'units','%');
x.subheading('');  
  x.series('Exchange Rate',exp(s.ls/100),'units','GEL/USD');
  x.series('',(s.ls-s.ls{-4}),'units','% (y-o-y)');
  x.series('Risk Premium',s.prem,'units','% of Exchange Rate');
x.subheading('');
  x.series('GDP',s.dot4_gdp,'units','% (y-o-y)');
x.subheading('');
  x.series('Interest Rate',s.rn,'units','% p.a.');

x.subheading('');
x.subheading('Monetary Conditions');
  x.series('Monetary Conditions',s.mci,'units','%');
  x.series('Real Interest Rate Gap ',s.rr_gap,'units','p.p.');
  x.series('Real Exchange Rate Gap',s.lz_gap,'units','%');  
    
x.subheading('');
x.subheading('Supply Side Assumptions');
  x.series('Potential Output',s.dot_gdp_eq,'units','% (q-o-q)');
  x.series('',(s.lgdp_eq-s.lgdp_eq{-4}),'units','% (y-o-y)');
  x.subheading('');
  x.series('Eq. Real Interest Rate',s.rr_eq,'units','%');
  x.subheading('');
  x.series('Eq. Real Exchange Rate',s.dot_z_eq,'units','% (q-o-q)');
  x.series('',(s.lz_eq-s.lz_eq{-4}),'units','% (y-o-y)');
    
x.pagebreak();   
    
x.table('Forecast Decomposition',TableOptions{:});

x.subheading('');
x.subheading('Inflation Decomposition');
  x.series('Inflation',s.dot_cpi,'units','%');
  x.series('Lag',s.b1*s.dot_cpi{-1},'units','p.p.');
  x.series('Expectations',(1-s.b1)*s.dot_cpi{+1},'units','p.p.');
  x.series('Real Marginal Cost',s.b2*s.rmc,'units','p.p.');
  x.series('RMC - Domestic',s.b2*s.b3*s.lgdp_gap,'units','p.p.');
  x.series('RMC - Imported',s.b2*(1-s.b3)*s.lz_gap,'units','p.p.');

x.subheading('');
x.subheading('Ouptut Gap Decomposition');
  x.series('Output Gap',s.lgdp_gap,'units','%');
  x.series('Lag',s.a1*s.lgdp_gap{-1},'units','p.p.');
  x.series('Monetary Conditions',-s.a2*s.mci,'units','p.p.');
  x.series('Real Interest Rate',-s.a2*s.a4*s.rr_gap,'units','p.p.');
  x.series('Real Exchange Rate',-s.a2*(1-s.a4)*(-s.lz_gap),'units','p.p.');
  x.series('Foreign Output Gap',s.a3*s.lx_gdp_gap,'units','p.p.');
  
x.subheading('');
x.subheading('Interest Rate Decomposition');
  x.series('Interest Rate',s.rn,'units','% p.a.');
  x.series('Lag',s.f1*s.rn{-1},'units','p.p.');
  x.series('Neutral Rate',(1-s.f1)*s.rn_neutral,'units','p.p.');
  x.series('Expected Inflation',(1-s.f1)*s.f2*(s.dot4_cpi{+4} - s.target{+4}),'units','p.p.');
  x.series('Output Gap',(1-s.f1)*s.f3*s.lgdp_gap,'units','p.p.');
  
x.pagebreak();   
    
x.table('Forecast - Foreign Variables',TableOptions{:});

x.subheading('');
x.subheading('USA');
  x.series('Inflation',s.dot_x_cpi,'units','% (q-o-q)');
  x.series('Interest Rate',s.x_rn,'units','% p.a.');
  x.series('Output Gap',s.lx_gdp_gap,'units','%');
        
x.publish('Forecast','display',false);
disp('Done!');